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Motivated by recent experimental and numerical studies of coherent structures in 
wall-bounded shear flows, we initiate a systematic exploration of the hierarchy of un- 
stable invariant solutions of the Navier-Stokes equations. We construct a dynamical, 
10^-dimensional state-space representation of plane Couette flow at Re ~ 400 in a small, 
periodic cell and offer a new method of visualizing invariant manifolds embedded in such 
high dimensions. We compute a new equilibrium solution of plane Couette flow and the 
leading eigenvalues and eigenfunctions of known equilibria at this Re and cell size. What 
emerges from global continuations of their unstable manifolds is a surprisingly elegant 
dynamical-systems visualization of moderate- i?e turbulence. The invariant manifolds tes- 
sellate the region of state space explored by transiently turbulent dynamics with a rigid 
web of continuous and discrete symmetry- induced heteroclinic connections. 



1. Introduction 



In a seminal paper, iHopi 



envisioned the function space of Navier-Stokes ve- 
locity fields as an infinite-dimensional state space, parameterized by viscosity, boundary 
conditions and external forces, in which each 3D fluid velocity field is represented as a 
single point. Laminar states correspond to equilibria that are globally stable for suffi- 
ciently large viscosity. As the viscosity decreases (Reynolds number increases), turbulence 
sets in, represented by chaotic state-space trajectories. Hopf 's observation that viscosity 
causes state-space volumes to contract under the action of dynamics led to his key con- 
jecture: that long-term, typically observed solutions of the Navier-Stokes equations lie on 
finite-dimensional manifolds embedded in the infinite-dimensional state space of allowed 
velocity fields. These manifolds, known today as 'incrtial manifolds,' arc well-studied in 
the mathematics of spatio-temporal PDEs. Their finite dimensionality for no n- vanishing 



viscos ity parameters has been rigorously established in certain settings by iFoias et 
and collaborators. 



Since Hopf's time, engineers and applied mathematicians have assembled a body of 
empirical evidence that moderately turbulent flow s exhibit organized, ii itrinsica l ly low - 
dim ensional behavior for a variety of conditions (see Holmes et all ( 19961 ). Panton ( 1997 ), 
and Robinsod |l99ll) for good overviews of this large body of work). The experiments 
of lKline et all (Il967f ). for example, revealed sp atially orga n ized s treaks in the turbulent 
boundary layer. The numerical simulations of Kim et al. ( 1987t ) opened access to the 
full 3D velocity field of channel flows and paved the way for mo r e deta iled studies of 
organization in wall-bounded flows. The work of [Hamilton et al. I (Il995l) began a very 
fruitful line of research; it identified from numerical simulations a remarkably well-defined, 
quasi-cyclic process among streamwise streaks and vortices (or 'rolls') in low-Reynolds 
number plane Couette flow. lWaleffel ( 19951 . 1997) further developed these ideas into a 'self- 
sustaining process theory' that explains the quasi-cyclic roll-streak behavior in terms 
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of the forced response of streaks to rolls, growth of streak instabilities, and nonlinear 
feedback from streak instabilities to rolls. 

The preponderance of recurrent, coherent states in wall-bounded shear flows suggests 
that their long-time dynamics lie on low-dimensional state-space attractors. This has 
motivate d a number of resea r chers to mode l such flows with low-dimensional dynamical 
systems. lAubrv et al. I (Il988l) : iHolmes et all I (|1996D used 'Proper Orthogonal Decomposi- 
tion' [POD] of experimental data and Galerkin projection of the Navier-Stokes equations 
to produce low-order models of coherent structures in boundary-layer turbulence. These 
models reproduce some qualitative features of the boundary layer, but the quantitative 
accuracy and the validity of simpli fying assump tions in their derivation are uncertain 
(|Zhou fc SirovichI ((1992 ): ,Sirovich fc Zhou 119941 : ICibsol] (|2002l )). POD models for plane 
Couette were developed by Smith et all ( 2005[ ) 

Another class of low-order mo dels of plane Coue tte flow derives fr om the 's elf-sustaining 
process theory' di sc ussed a b ove (lOauchot fc Viouiard (,2000, ) ; .Moehlis et al. (|2004 i2005il ; 
Mannevilld (l2004l ): ISkufcal tOO^ ). These models use analytic basis functions explicitly 



designed to represent the streaks, vortices, and instabilities of the self-sustaining process, 
compared to numerical basis functions of the POD, which represent statistical features 
of the flow. They improve on the POD models by capturing the linear stability of the 
laminar flow and saddle- node bifurcations of non-trivial 3D equih bria consisti ng of rolls, 
streak s, and streak undulations. The work of lSkufca et al. ( 20061 ). based on a ISchmiegel 
( 19991) 9-variable model, offers an elegant dynamical systems picture, with the stable 
manifold of a periodic orbit defining the basin boundary that separates the turbulent 
and laminar attractors at Re < 402 and the stable set of a higher-dimensional chaotic 
object defining the boundary at higher Re. However, these models share with POD models 
a sensitive dependence on modeling assumptions and uncertain quantitative relations to 
true Navier-Stokes flows. A systematic study of the convergence of POD/Galcrkin models 
of plane Couette flow to fully-resolved simulations indicates that dimensions typical in 
the literature (10-10^) are orders of magnitude too low for either short-term quantitative 
prediction or reproduction of long-term statistics ([Gibson (2003 )). 

The lack of quantitative success in low-dimensional modeling motivates yet another 
approach: the calculation of exact invariant solutions of the fully-resolved Navier-Stokes 
equations. The idea here is to bypass low-dimensional modeling and to tre at fully- 



resolv ed CFD algorithms directly as very high-dimensional dynamical systems. iNagatal 
( 199(1 ) computed a 'lower-branch' and 'upper-branch' pair of nontrivial equilibria solu- 
tions to plane Couette flow by continuation and bifurcation from a wavy vortex solution 
of Tayl or-Couette flow. Starting w ith physical insights from the self-sustaining process 
theory, IWaleffd (Il998l . l200lL boOSi) generated, ab initio, families of exact 3D equilibria 
and traveling waves of Navier-Stokes in plane Couette and Poiseuille flows for a variety 
of boundary conditions and Re numbers, using a lO'^-dimensional Newton search and 
continuation from non-equilibrium states that approximately balanced the mechanisms 
hig hlighted by the self-sus tained process. As noted in IWaleffe C2Q03i), these solutions, 
and I Clever fc Busse (1992)'s equilibria of plane Couette flow with Rayleigh-Benard con- 
vection, are ho motopic to the Nagata eq uilib ria under smooth t ransfo rmations in the 
flow conditions. Faisst fc Eckhardt ( 2003 ) and Wedin fc Kerswell (2004) carried the idea 
of a self-sustaining process over to pipe flow and applied Waleffe's continuation strat- 
egy to derive families of trave ling-wav e solu tions for pipes. Traveling waves for plane 
Couette flow were computed bv lNagatal(|l997t ) using a contin uation method. Late r, trav- 
eling waves for pressure-driven channel flow were obtained bv lltano fc Toh (2001) with a 
shooting met hod. The first short - perio d unstable periodic solution of Navier-Stokes were 
computed bv lKawahara fc Kidal (|200lt l. Recentlv. Iviswanathi (|20076l ) has computed rel- 
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ative periodic orbits (orbits which repeat themselves with a translation) and further 
periodic orbits of plane Couette flow that exhibit break-up and reformation of roll-streak 
structures. 

The exact solutions described above turn out to be remarkably similar in appearance 
to coherent structures observed in DNS and experiment. iWalcffc (,2001ll coined the term 
'exact coherent structures' to emphasize this connection. The upper-branch solution, for 
example, captures many statistical features of turbulent plane Couette flow and appears 
remarkably similar to the roll-streak structures obse rved in direct numerical simulations 
(compare ?^g\lve^h) to fignve^c)). IWaleffd (|2003f ) showed that the upper and lower- 
branch equilibria appear at lowest Reynolds n umber with streak spacing of 100^ wall 
un its, an excellent match t o tha t observed in Kline et al\ (|l967l ). The periodic orbits 
of Kawahara fc Kida ( 2001 ) and Viswanath ( 200761 ) appear to be embedded in plane 
Couette flow's natural ergodic measure, and most of them capture basic statistics more 
closely than the equilibria. In pipe flow, high speed streaks that match the traveling-wave 
solutions in cross-section h ave been observed in beautiful e x periments using stereoscopic 
particle imag e vclocimctry (iHof et al\ (12004 ) : IBuss^ (|2004|) : Barenghil (|2004[) ). Addition- 
ally, there is preliminary evidence that the instabilities of these exact solutions play 
important dynamic roles. The relevance of stead y solut ions t o sustained tu r bulenc e and 
transition to turbulence is discussed in IWaleffd (|2003l ) and Ijimenez et tZI (|2005[ ). The 
stable manifold of the lower-branch solution is con jectured to control the b asin bound- 
ary be t ween the tu r bulent a nd laminar attractors (Waleffe fc Wand ( 20051 ): Wang et al 



(|2007f ): IViswanathI (|2007q[ )). iKerswell fc Tuttvl (|2007[ l 's numerical simulations have es- 
tablished that lower-branch traveling waves as act as similar boundaries in pipe flow, and 
that turbulent fields make occasional visits to the neighborhoods of traveling waves. 

Together, these results form a new way of thinking about coherent structures and tur- 
bulence: (a) that coherent structures arc the physical images of the flow's least unstable 
invariant solutions, (b) that turbulent dynamics consists as a series of transitions be- 
tween these states, and (c) that intrinsic low-dim ensionali t y in t urbulence results from 
the low number of unstable modes for each state ( Waleffe ( 20021 )). The long-term goals 
of this research program are to develop this vision into quantitative, predictive descrip- 
tion of moderate- i?e turbulence, and to use this description to control flows and explain 
their statistics. Much of this has already been accomplished in the simpler context o f the 



Kuramoto-Sivashinsky equation (jChristiansen et al\ (|1997l ) ; ICvitanovic et al\ ()200 



In this paper, we take a few steps towards realizing these goals in the case of plane 
Couette flow. In §[5] we review the physical characteristics and symmetries of plane Cou- 
ette flow. §[3] discusses the computation of invariant solutions and their eigenvalues and 
presents (a) a new equilibrium solution of plane Couette and (b) the linear stability 
analysis of this and the lower and upper-branch equilibria. These computations set the 
stage for the main advance reported in this paper, visualization and exploration of the 
state space of moderate- i?e plane Couette flow, undertaken in §|4l The combination of 
equilibrium solutions, linear stability analysis, and state-space portraiture reveals pre- 
viously unseen dynamical connections amongst the known invariant solutions of plane 
Couette flow. Particularly beautiful and unexpected arc the discrete symmetry enforced 
interrelations between unstable manifolds manifest in figure [S]- figure [9l 



2. Plane Couette flow 

Plane Couette flow is comprised of an incompressible viscous fluid confined between 
two infinite parallel plates moving in-plane at constant velocities. We take the length 
scale L to be half the distance between the walls and the velocity scale U to be half the 




(a) (b) 

Figure 1. Snapshots of plane Couette turbulence at .Re — 400. Velocity fields {u,v,w) 
in a periodic cell of size n = [77T/4, 2, 67t/5] (|Hamilton et al.\ (fl995h ') are shown with arrows for 
in-plane velocity and a colormap for the streamwise velocity component u: red/blue indicates 
u = ±1; green, u = 0. The upper wall at j/ = 1 and the upper half of the fluid is cut away 
to show the velocity in the y — midplane. The two snapshots shown are different instants 
from a simulation initiated with a random pertubation, selected to show (a) minimum and (b) 
maximum organization in the turbulent field. In particular, (b) resembles the upper-branch 
equilibrium shown in figure [Jl'cJ. 



relative wall velocity. After nondimensionalization and absorption of fluid density into 
the pressure field, the Navier-Stokes equations take the form 

d\x 1 

— -fu-Vu=-VpH V^u, V-u = 0, 

ot Re 

where the Reynolds number is defined as Re = UL/iy and v is the kinematic viscosity 
of the fluid. The plates move at speed ±1 along the 'streamwise' x-axis, the direction 
normal to the plates is the 'wall-normal' j/-axis, y G [—1,-1-1], and the in-plane z-axis, 
normal to the plate velocity, is referred to as 'spanwise.' The x,y,z unit vectors are 
x,y, z. (We use boldface to indicate vectors in three spatial dimensions.) The velocity 
field u has streamwise, wall-normal, and spanwise components u = [u, v, w]; the velocity 
at point X and time t is u(x, i) = [u,v,w]{x,y, z,t). The no-slip boundary conditions at 
the walls are u(x, ±1, z) = [0, ±1, 0]. Numerical computations replace the infinite x and 
z directions with a periodic cell of lengths and , or equivalently, the fundamental 
wavenumbers a = 2n/Lj. and 7 = 2n/Lz. We denote the periodic domain of the cell by 
f2 = [0, X [—1,1] X [0, Lz] or simply VI = [L^,, 2, L^]. We assume that the spatial mean 
of the pressure gradient is zero, i.e. there is no pressure drop across the cell in x or z. 

Replacing u with u -I- y x recasts Navier-Stokes in terms of the difference of the velocity 
from laminar flow: 

— -^y— -Hwx + u- Vu = -Vp+— V^u, V-u==0. 2.1) 

ot ox Re 

The difference u satisfies Dirichlet conditions at the walls, u(a;,±l,2) = 0. Henceforth 
we refer to the difference u as 'velocity' and u-|-yx as 'total velocity,' and we take (|2.ip 
as the Navier-Stokes equations for plane Couette flow. 

Plane Couette flow is the simplest of all shear flows, and it is here that roll-streak 
structures take their simplest form. For moderate values of i?e, the rolls span the full 
distance between the walls, whereas in channel and boundary-layer flows such structures 
are bounded by a wall on one side and open flow on the other. Figure [1] shows two typical 
velocity fields from a simulation in the 'HKW cell ^ = [77t/4, 2, 67t/5] at Re = 400. 
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The numerical simulations of iHamilton et al. ( 1995 ) indicate that this is roughly the 
smallest cell and Reynolds number that sustains turbulence for long time scales. Roll- 
streak structures arc evident, particularly in figure [T]f'6j. The rolls circulate high-speed 
fluid towards the walls and low-speed flow away; the resulting streaks of high-speed fluid 
near the walls dramatically increase drag compared to laminar flow. For example, the 
power input needed to maintain constant wall velocity in plane Couette flow increases by 
a factor of three if the flow goes turbulent (sec figure [21(e)). The practical importance of 
roll-streak dynamics derives from their role in momentum transfer and turbulent energy 
production and their generic occurrence in wall-bounded shear flows. 

Except for figure [1] and parts of figure [2l the re s ults in this paper arc for Re = 400 
and n = [27t/1.14,2,47t/5], first studied in Twaleffel (|2002l ). This cell matches the HKW 
cell [77t/4, 2, 67t/5] closely in x (7/4 w 2/1.14). The z length scale Lz = 47t/5 was chosen 
as a compromise between Lz = 67t/5 of the HKW cell (which sustains turbulence for 
long time scales but has equilibria only with doubled period in z) and its first harmonic 
Lz = 37t/5 (which has equilibria at the fundamental harmonic, but tends to decay to 
laminar flow). Simulations for these parameters tend to decay to the laminar state in 
within several hundred nondimensionalized time units L/U , but the transient dynamics 
serves well to illustrate our invariant manifolds construct ion. Whether a given cell size 
sust ains turbulence in definitely is a subtle dynamical issue: ISchmiegel fc Eckhardd (|l997h 
and lSchmiegel (|l999l ) observe only chaotic transients in their studies. 



2.1. Energy transfer rates 

The kinetic energy density E, the bulk dissipation rate D, and the power input / of total 
velocity field of plane Couette flow are given by 



E{1^^^ jjt^]^\vi + yk\^ 

D{t) - :^ /dx|V X (u + yi) 
V Jn 



^(*) = l + 7^ / dxdz 



2A 



di 



dy y=i dy 



(2.2) 
(2.3) 
(2.4) 



where V = iL^Lz and A = L^Lz- The normalizations are chosen so that D = I = 1 
for laminar flow and E = I — D. Yigxvce^e) shows / vs. D for a turbulent trajectory 
and several invariant solutions. Equilibria and relative equilibria must fall on the line 
I = D where these two quantities are in balance. The energy input and dissipation rates 
must also balance in averages over a single period of a periodic orbit or relative periodic 
orbit p. Dp = 1/Tp J^''dt D{t) ~ 1/Tp ^^^dt I{t) ~ Ip, as well as for long-term averages. 

We note that the turbulent trajectory in figure O stays clear of the much lower dis- 
sipation rates Dlb = 1-429 and I^nb = 1-454 of the equihbrium solutions Ulb and Unb 
(see ? l3.2p . so these equilibria are far from the turbulent a ttractor . The e nergy, the dis- 
sipation rate, and mean and RMS velocity p rofiles of the iNagatal (|l990f) upper-branch 
equilibrium and the iKawahara fc Kidal (j20011 ) periodic orbit are all numerically close to 
the long-term turbulent averages. This suggests that the solutions play an important role 
in turbulent dynamics, but turbulent statistics do not simply follow from the properties 
of one or two solutions. On the contrary, periodic orbit theory shows that the statistics 
of dynamical systems are given by sums over hie rarchies of periodi c orbit s, with weights 
determined by the orbits' lengths and stabilities ( Cvitanovic et al. ( 2007t l). 
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Figure 2. (a,b) Spatial-mean and spatial-RMS velocity profiles of the iNagatal (| 19901 ) up- 
per-branch equilibrium (symbols) in S7 = [27t/1.14, 2, 47t/5] against temporal mean and RMS 
profiles for sustained turbulent DNS data in [77t /4, 2, 67t/5] (lines). Re = 400 for both. (c,d) 
Temporal mean and RMS velocity profiles for the iKawahara fc Kidal (|200ll ') Re = 400 periodic 
solution (symbols), against the same turbulent data shown in (a,c). In (a,c), mean values of u 
are marked with o; in (b,d) RMS values of u, v, w are marked with o, A, □. (e) I versus D for a 
typical trajectory of sustained turbulence in the [77t/4, 2, 67t/5] cell (lines), mean dissipation rate 
over IKawahara fc Kidal (|2001l ') periodic orbit D = 2.91 (o), mean turbulent dissipation D = 2.89 
(A), and upper-branch equilibrium Dub = 2.89 (□). The laminar equilibrium dissipation rate 
is -D = 1 (not shown). 



2.2. Symmetries 

Plane Coucttc flow is invariant under two discrete symmetries ai , a2 and a continuous 
two-parameter group of translations t(£x,£z)'- 

Ui [u, V, w] (x, y, z) = [u, V, ~w\ {x, y, ~z) 

(72 [u,v,w\{x,y,z) = [-u,~v,'w\{-x,-y,z) (2.5) 
t(4,4)['«,w, w](a;,y,z) = [u, w, w;](x -I- 4, y, z -I- 4) • 

The Navier-Stokes equations and boundary conditions are invariant under any symmetry 
s in t he group gener ated by these symmetries: d{su)/dt — s{du/dt). 

The lNagatal (|l990[ ) lower and upper-branch equilibria Ulb and Uub are invariant under 
action of the subgroup S* = {1, si, S2, S3}, where si = t{Lx/2, 0) (Ti, S2 = t[Lx/2^ Lz/'i-) C2, 
and S3 = S1S2. That is, s Ulb = and s Uub = Uub for s G S. The si and S2 symmetries 
are referred to as the 'shift-reflect' and 'shift-rotate' symmetries. The group actions on 
velocity fields u are given by 



si [m, V, w]{x, y, z) = [u,v, -w]{x + Lx/2, y, -z) 

52 [u,v,w]{x,y,z) = [-u,-v,w]{-x + La;/2, ~-y, z + 1^/2) 

53 [u, V, w] [x, y, z) = -w, -w] (-x, -y, ~z + 1^/2) . 



(2.6) 



We denote the space of velocity fields that satisfy the kinematic conditions of plane 
Couette flow by 

U = {u I V • u = 0, u(x, ±1, z) = 0, u(a;, z) = u{x + L^, y, z) = u{x, y,z + L^)} ■ 

(2.7) 

and the 5- invariant subspace ( Golubitskv fc Stewart ( 20021 ) ) of U by 

Us = {u G U I sju = u , Sj G 5} (2.8) 
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Us is a flow-invariant subspace of U since S symmetry is preserved by evolution under 
the Navier-Stokes equations. 

A second important subgroup is the group of half-cell translations T = {1, r^, t^z}, 
where = t{Lx/2^0), Tz = r(0, and Txz = t^Tz- In general, the continuous 

translation r^ix^iz) maps each state u into a 2-torus of dynamically equivalent states, 
and the group {1, cri, (T2, crifT2} maps these into four dynamically equivalent 2-tori. For 
u S U5, the four tori coincide, and the torus intersects U5 at the four points ru, r G T. 
(Since elements of T commute with those of 5, u g U5 implies tu e U5 for r G T.) For 
example, the upper-branch equilibrium Uub appears within U5 in four distinct half-cell 
translations, namely Uubj 'HeUub: TzUu-b, and T^zUuB- 



3. Invariant solutions of plane Couette 

Let F(u) represent the Navier-Stokes equation (|2.ip for u e U (|2.7p and f* its time-t 
forward map 

F(u), f*(u) = u+ rdrF(u). (3.1) 



F(u) admits of invariant solutions of the following types: 
F(ueq) = equilibrium or steady state Ubq 

F(utw) = • Vutw relative equilibrium or traveling wave u^w, velocity c 
f (up) = Up periodic orbit p of period Tp 

f'^p(up) = Tp Up relative periodic orbit, period Tp, shift Tp = T{£x,iz) ■ (3.2) 

Relative equilibria and relative periodic orbit solutions are allowed due to the continuous 
translation symmetry T{£a:,£z)', for traveling waves, boundary conditions require c-y = 0. 
We expect to see many more relative periodic orbits than periodic orbits because a 
trajectory that starts on and returns to a given torus is unlikely to intersect it at the 
initial point, unless forced to do so by a discrete symmetry. This indeed is the case for 
other PDEs with co ntinuous symmetries, such as the complex Ginzburg-Landau equatio n 
( Lopez et all (|2006l )) and the Kuramoto-Sivashinsky equation (jCvitanovic et 



Restriction to the S'-invariant subspace Us defined in (|2.8p eliminates relative equilibria 
and relative periodic orbits. In what follows we focus mostly on dynamics within Us. 

3.1. Finite representation 

Computing the exact solutions and stability modes of plane Couette flow requires a 
finite but fully-resolved discretization of the constrained partial differential and integral 
equations represented by (|3.ip and (|3.2p . We investigated two approaches to discrete 
representation. In the first approach the vector u G was formed by breaking the 
complex spectral expansion coefficients of a CFD algorithm into real and imaginary 
parts and then selecting from these a set of linearly independent real- valued coefficients. 
Our CFD algorithm, chann elf low . org, is based on the velocity-pressure algorithm of 
Kleiser fc SchumannI (|l980( )) with expansions 



u(x, <) = X] E E E "^-fc'™ ^^(y) 6^™^^"^/^^+'=^/^^) , (3.3) 

j = - J k=-K e=0 m=l 

where the Te are Chebyshev polynomials and (xi,X2,X3) = (x, y,z) unit vectors. The 
algorithm employs a Chebyshev tau method and tau correction for enforcement of incom- 
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pressibility and boundary conditions, third-order semi- implicit backwards-difFerentiation 
time-stepping, dealiasing in the x, z transforms, and a variety of methods for calculating 
the nonlinear term. The expansion (|3.3p retains a number of linearly dependent terms, 
due to complex symmetries and the run-time enforcement of the incompressibility and 
boundary conditions. Intimate knowledge of the CFD algorithm and careful accounting 
is required to determine the precise value of the dimension of the linearly independent set 
and a self-consistent method of converting back and forth between the state-space vector 
u and the expansion coefficients of u. For our CFD algorithm and an iV^ x Ny x N^^ grid, 
d is slightly less than 2Nx{Ny — 2)Nz. The accounting for velocity- vorticity algorithms 
is somewhat simpler since incompressibility is elimin ated a t the outset. For further tech- 
nical details, please refer to Iviswanath (,20076. ). .Halcrowl ( 20081 ). and documentation at 
channelf low . org, or contact the authors. 

In the second approach, we explicitly constructed a set of orthonormal, divergence-free, 
no-slip basis functions ^n(x) and formed the state-space vector u from the coefficients 
Un of the expansion u(x) = 'Yli'L=i^^^^i^)- This approach produces a mathematically 
simpler representation, in that (1) all constraints are subsumed into the basis and elimi- 
nated from further consideration, (2) an explicit second-order ODE of form ii = F{u) can 
be derived through Galerkin projection of the Navier-Stokes equation onto the basis set, 
and (3) with proper normalization of the basis functions, the norm of the state-space 
vector u e is the same as the energy norm of the velocity field u. The downside is 
that the formulation of the basis set is complicated, and it requires extra computation 
for orthogonalization and transforms between state-space vectors and the CFD represen- 
tation. We found no practical advantages to the orthonormal basis. The results reported 
here were computed using CFD expansion coefficients for the state-space vector u. 

The choice of discretization u GMf^ and CFD alg orithm implicity defines a d-dimensional 
dynamical system it — F{u). IViswanathI (20076) showed invariant solutions and linear 



stability of F can be computed efficiently with Krylov subspace methods and numer- 
ical evaluation of the finite-time map : u{t) ^ u{t + T) with the CFD algorithm. 
Equilibria may computed as solutions of f^{u) ~ u = for fixed t; and periodic orbits 
as solutions of the same equation with varying t. Viswanath's algorithm for computing 
these solutions involves a novel combination of Newton descent, GMRES solution of the 
Newton equations, and 'trust-region' limitation to the magnitude of the Newton steps. 
The results reported in this paper, however, used straight Newton-GMRES search, with 
no trust region modification. We will often discuss equilibria and linear stability in terms 
of the fiow F, with the understanding that the computations are performed using the 
finite-time map f^. 

3.2. Equilibria 

The starting point s of our exploration of plane Couette state space arc the Nagata ( 1990[ ) 
and lWaleffd ( 20031 ) Uub and Ulb equilibrium solutions of plane Couette flow for Re = 400 
and [Lx, Ly,Lz\ = [27t/1.14, 2, 47t/5], provided in numerical form by Walcffe. These fields 
employed an elliptical truncation of spectral coefficients (requiring P / + /K^ + 
/LF' < 1 for expansions of the form (|3.3p ) on a 32 x 34 x 32 grid. We use Viswanath's 
Newton-GMRES algorithm to increase the resolution to a rectangular truncation < 
-^,1^1 < K,l < L) on a 32 X 35 X 32 grid. The dimensionality of this finite representation 
is c? = 61 506. These recomputed solutions satisfy (|3.2[) 



F{u)^0, u = f{u) (3.4) 

in discrete form. At this spectral resolution, the Newton-GMRES search can reduce the 
residual of the discrete equilibrium equations to 10~^*, but the truncated coefficients are 
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(a) (b) (c) 

Figure 3. Equilibrium solutions of plane Couette flow: (a) Ulb, the lower-branch 
equilibrium, (b) Unb, the 'newbie' equilibrium, and (c) Uub, the upper-branch equilibrium. 
[L^c, Ly, L-\ — [27t/1.14, 2, 47t/5] and Re = 400. The plotting conventions are the same as in 
figure [1] 



of the order 10~^, and the residual increases to the same level when the given solution is 
integrated at higher resolution. The accuracy of the equilibrium solutions is thus roughly 
single-precision. 

Figure [3] shows the Ulb and Uub equihbria as velocity fields, alo ng with a third equili b- 
rium, Unb, that was discovered in the course of this investigation ( Halcrow et al. ( 20081 )). 
We refer to the new equilibrium as u^b, pronounced 'ncwbic,' in keeping with the nomen- 
clature for other equilibria. This equilibrium was found by initiating Newton-GMRES 
searches for zeroes of the equilibrium equations from points within the unstable mani- 
folds of Ulb and Uub (sce §|1]) . A portion of the one-dimensional unstable manifold of Ulb , 
shown in figure O appears to be strongly influenced by a complex unstable eigenvalue of 
Unb- Initial guesses along this portion of the Ulb unstable manifold converge rapidly to 
Unb J to several digits of accuracy in a few Newton steps. 



3.3. Linear stability of equilibria 

Dynamics in the neighborhood of an equilibrium solution are governed by the linear 
stability matrix 

[DFU^ = ^ (3.5) 

Figure [4] shows the lea ding (most unstabl e) eigenvalues of Ulb, u^b, and Uub, computed 
with Arnoldi iteration ( ViswanathI ( 200761 )). Figure Ulf'aJ shows all computed eigenvalues; 
figure IH'S j shows those within the 5- invariant subspace Ug. The eigcnfunctions v of DF 
at the equilibria are either symmetric, sv — v, or antisymmetric, sv = — v, as all s S S' 
are idempotent, = 1. Thus, in general, the dynamics carries small perturbations of 
these equilibria into the full space U. The Ulb equilibrium has a single unstable eigenvalue 
(jWang et (120071) 1. Within Us, the Uub equilibrium has a single unstable complex pair, 
and Unb has one unstable real eigenvalue and one unstable complex pair. 

The Arnoldi eigenvalues are accurate to 10~^, as determined by repeated calculations 
with different random initial vectors, and comparison of Arnoldi computations to an- 
alytically known eigenvalues of the laminar equilibrium. This level of accuracy results 
from our use of off-center finitc-diffcrcncing to estimate differentials of the flow in the 
Arnoldi iteration: Df*\uV = {f{u + ev) - f\u))/e + 0(e), with e = lO"'^. Tables of 
numerical eigenvalues and their symmetries are given in §E] (table [T] and HD and at 
channelf low . org. 
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Figure 4. Leading Ulb, Unb, Uub eigenvalues in (a) the full space U (b) the S- invariant 
subspace Ug. Numerical values are listed in tables [1] and [2] 



3.4. Linearized evolution 

Let A,?;eq be an eigenvalue, eigenvector solution of DFIu^qV = Xv at the equilibrium 
WeqI3 Then the linearized state-space dynamics v = DFIu^qV about u^q has solution 
v{t) = e^'^VEQ, and the initial condition u(0) = Ueq + ^Veq with e \veq\ ^ 1 evolves as 

= UEQ + ei'EQe^* + 0(e2). (3.6) 

The linearized evolution of the velocity field u(x, t) can then be derived by reconstructing 
the velocity fields from the corresponding state-space vectors, as discussed in § 13.11 Small 
perturbations about Ueq along the cigenfunction Veq evolve as 

u(x,i) =UEQ(x) + evEQ(x)e^* + 0(e2). (3.7) 

Complex eigenvalues and eigenvectors must be recast in real-valued form prior to 
conversion to velocity fields, since each element of the state-space vector v is the real or 
imaginary part of a complex- valued spectral coefficient in a CFD expansion such as (|3.3[) . 
Let Aeq'""''^'' = ± iw^") be a complex eigenvalue pair and wiq'"^"'^'' = Vr^^ ± w^"' the 
corresponding complex eigenvectors. Then (dropping superscripts) the initial condition 
u(0) = Ueq + evr cvolvcs as a real-valued spiral 

u{t) = Ueq + e {vr cosLut — Vi siuwt) 6^* + O(e^) . (3.8) 

Real- valued fields v,, and v; can be reconstructed from the real- valued vectors and Vi , 
and the real- valued initial velocity field u(0) = Ueq + e evolves as 

u(i) = Ueq + e (vr coswt - v; sinwi)e^* + O(e^) . (3.9) 

3.5. Unstable manifolds 

Let (Weq) denote the stable (unstable) manifold of equilibrium Ueq. For each real- 
valued unstable eigenvalue A'"\ wc shall refer to the orbit of an infinitesimal perturbation 
of Ueq along the corresponding cigenfunction Veq^ as Weq'\ This part of the Ueq un- 
stable manifold is 1-dimensional and can be computed by DNS integration of the initial 
conditions Ueq ± cVeq^ where e ^ 1. 

For an unstable complex pair A^"'"^^-' of equilibrium Ueq, let II^eq"'"'''^^ denote the 
orbit of a circle of infinitesimal radius in the plane about Ueq spanned by vi"'' , v^^"-* . 

f We indicate particular invariant solutions with subscripts, such as ulb or Ulb for the 
lower-branch equilibrium solution. The nth eigenvalue is A'"' , n = 1, 2, . . ., in order of decreasing 
real part. Whenever the context allows it, we shall omit the eigenvalue and/or solution labels. 
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This part of the Ueq unstable manifold is 2-dinicnsional; its shape can be traced out by 
computing a set of trajectories with initial conditions Ueq + e(vr"' cosO + vl"-* sin 6*) for 
a set of values of 9. In practice, one obtains a more uniform distribution of trajectories 
by setting initial conditions along the line Ueq + evf"', for a set of values of e. 

The global unstable manifolds W^q'^ and M^eq"'"^^'' arc invariant sets that preserve 
the symmetries shared by the equilibrium and the eigenvectors from which they are 
generated. The S'-invariant subspace portions of the unstable manifolds of Ulb, Uubi and 
Unb have dimensionality of 1, 2, and 3, respectively, see figure In what follows, we 
will focus on these low-dimension unstable manifolds confined to the Us subspace. 



4. The geometry of plane Couette state space 

We now turn to the main theme of this paper: exact state-space portraiture of plane 
Couette flow dynamics. The state-space portraits are dynamically intrinsic, since the 
projections are defined in terms of solutions of the equations of motion, and represen- 
tation independent, since the projection operation (the inner product (j4.ip ) is indepen- 
dent of the numerical representation. The method is by no means restricted to plane 
Couette flow or our choice of state-space representation or CFD algorithm. It can be 
appli ed to any high-dimensiq nal dissipative flow, for example the Kuramoto-Sivashinsky 
flow ( Cvitanovic et al. ( 2008f l). Production of state-space portraits requires numerically 



computed physical states (such as equilibrium solutions and their linear stability eigen- 
functions), an algorithm for integrating dynamics, and a method of computing the inner 
product between states over the physical domain. 

4.1. Peering into oo- dimensional state spaces 

Numerical methods have advanced to the point where it is possible to compute highly 
accurate unstable exact coherent states in low-Reynolds shear flows. How is one to visu- 
alize them? Even though fully-resolved solutions of Navier-Stokes equation are embedded 
in 10^ or higher dimensional state spaces, there are few unstable eigendirections for Re 
close to the onset of turbulence. The associated asymptotic strange attractors / repellers 
might thus be amenable to dynamical systems visualizations, such as trajectory projec- 
tions, Poincare sections, state-space partitions and symbolic dynamics description. 

In this section, we show that revealing, representation-independent state-space por- 
traits can be deflncd in terms of invariants of the dynamical system. The idea is to 
choose as a basis set states of the fluid with characteristics of recurrent coherent struc- 
tures, and to project the evolving fluid state u{t) onto this basis with the energy norm 
(|2.4p inner product 

(u,v) = / dx u- V, llujp = (u,u) . (4.1) 

That is, we form orthonormal basis functions {ei, 62, . . . ,e„} from a set of linearly 
independent fluid states and produce a state-space trajectory 

a(t) (ai,a2,--- ,a„,- ••)(*), a„(t) = (u(t), e„) (4.2) 

in the {e„} coordinate frame by (|4.ip . The projection can be viewed in any of the 2d 
planes {em,e„} or in 3d perspective views {e£,em,e„}. The resulting portraits depend 
on the fluid states involved and not on the choice of numerical representation. Orthonor- 
mality of the basis set is not strictly necessary, but with it, distances are directly related 
to (|2.4|) . the energy norm of u. 
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The low-dimensional projections presented in this section are closely related to other 
finite approaches to Navier-Stokes in both technical methods and purpose. For example, 
on a technical level, the projections in this section differ from the finite discretizations 
discussed in ii l3.1l only by degree of dimensionality. If the dimension n of the state- 
space representation (j4.2[) were taken to the dimension d of the fuUy-resolvcd numerical 
discretization, the two discretizations would be related by a simple linear transforma- 
tion. We emphasize the differences between the two through notation: u for the high- 
dimensional vector of coefficients of a fully-resolved numerical discretization, and a for 
the low-dimensional coordinates of a state-space portrait. 

The proj ection methods her e are quite similar in spirit to the low-dimensional projec- 
tions of the lAubrv et al\ (|l988l ) POD modehng approach, in that they aim to capture key 
features and dynamics of the system in just a few dimensions. Indeed, our use of the 
inner product, orthonormal basis functions, and the very idea of constructing a basis from 
characteristic states derive directly from POD modeling. But the methods presented here 
depart from the POD in two key points: (1) We construct basis sets from exact solutions 
of the full-resolved dynamics rather than from the empirical eigenfunctions of the POD. 
Exact solutions and their linear stability modes (a) characterize coherent fluid states 
precisely, compared to the truncated expansions of the POD, (b) allow for different basis 
sets and projections for different purposes and different regions of state space, and (c) 
are not limited to Fourier modes and 0(2) symmetry in homogeneous directions. (2) 
We deploy low-dimensional visualization without any low-dimensional modeling. The dy- 
namics are computed with fully-resolved direct numerical simulations and projected onto 
basis sets to produce low-dimensional state-space portraits, tailored to specific purposes 
and specific regions of state space. The portraits reveal dynamical information visually, 
providing insight to dynamics that can guide further analysis. Specifically, we do not 
suggest that any of our low-dimensional projections is suited to a global projection of 
the state-space dynamics into a low-dimensional ODE model. 



4.2. A global basis spanned by discrete translations of u,jg 

There is an infinity of possible basis sets, but two choices appear particularly natural: 
(a) global basis sets, determined by a set of dynamically important and distinct states, 
or (b) local basis sets, defined in terms of a given equilibrium Ueq and its linear stability 
eigenfunctions Veq^ . An example of a local coordinate system based on eigenfunctions of 
the equilibrium is presented in § 14.41 an example of a global basis is defined here and 
used to construct state-space portraits in § 14.31 

The projection for a global state-space portrait should emphasize important global 
features of the flow. For example, for a system with three distinct equilibria, a good flrst 
guess for a plane of projection would be the plane containing the three equilibria. The 
system under study has three distinct equilibria Uub, Ulb, and Unb, each appearing in 
four spatial phases, plus the laminar equilibrium at the origin. We have found that for the 
S-invariant subspace Us the irreducible representations of the half-cell translations group 
T (§ 12. 2p provide natural linear combinations of a given equilibrium and its translations. 

For example, a set of orthonormal basis functions based on Uub and its half-cell trans- 
lated siblings can be generated by the four irreducible representations of the D2 dihedral 
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group T = {l,Tx,Tz,Txz} (see §[22]): 
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(4.3) 



where c„ is a normalization constant determined by j|e„|| = 1. The last 3 columns indicate 
the symmetry of each basis function under half-cell translations; e.g. S in the column 
implies that r^en = e„ and an A that rxe„ = — 6„. As the 'velocity' u in the Navier- 
Stokes equation (|2.ip for plane Couette flow is the difference from laminar flow, the origin 
in state-space portraits corresponds to the laminar equilibrium Ulm • This origin is shared 
by all symmetry-invariant subspaces, as Ulm = is invariant under all symmetries of the 
flow. Note, however, that the basis functions e„ are not themselves invariant solutions 
of Navier-Stokes; rather, they form an orthogonal coordinate system that spans the four 
translations of Uub within the S*- invariant subspace Vs- 

The evolution of a state u € Us is represented in this projection by the trajectory 
a{t) = (ai, 02, 03, a4)(t) with a„(i) = (u(i),e„). As discussed in ii l4.11 this is a low- 
dimensional projection intended for visualization. The dimensionality is lower than the 
full state space, so trajectories can appear to cross in such projections. We emphasize 
again that this is one of many possible projections that can be constructed from linear 
combinations of exact solutions, their spatial translations, and their eigenfunctions. An 
example of a more complex basis construction is given in ii l4.4l 



4.3. A global stroll through plane Couette state space 

With this road map in hand, let us take a stroll through the state space of a transiently 
turbulent plane Couette flow. Like many dynamical narratives, this will be a long walk 
through unfamiliar landscape with many landmarks of local interest. We undertake the 
tour for several reasons. The main message is that now such a promenade is possible 
even in 10^ dimensions. But a detailed road map is a necessary prerequisite for solving at 
least three outstanding problems: (a) uncovering the interrelations between (in principle 
infinite number of) invariant solutions, such as those of figure El (b) a partition of state 
space is a needed for a systematic exploration of dynamical invariant structures such as 
relative periodic orbits, and (c) explicit linear stability eigenvectors and their unstable- 
manifold continuations will be needed to control and chaperone a given fluid state to a 
desired target state. 

Our first example of a global state-space portrait of plane Couette flow is flgure [5l 
Here trajectories in the unstable manifolds of Ulb, u^b; u^b and several of their half- 
cell translations arc projected onto {ei,e2} plane defined by (|4.3p . Both ei and 62 are 
symmetric in t^, so points related by half-cell translations in x (such as Ulb and t^Ulb) 
map to the same point in this projection. The basis function 62 is antisymmetric in , so 
half-cell translations in z appear symmetrically opposite along 02- u^b and its unstable 
manifold are shown only in a single z translation, since the reversed orientation of the 
unstable spiral of TjjUnb only clutters the picture. 

The Ulb unstable manifold W^^ forms the backbone of the dynamics shown in 
Figure [5l The Ulb equilibrium has a single real -valued unstable eig envalue, as shown in 
figure [4] for Re = 400 and for much higher Re in Wang et al. ( 2007[ ). W^b therefore Id 



and can be computed in its entirety as discussed in §[331 Figure [5l shows the projection 
of onto {61,62} plotted with thick blue lines. Both branches of decay to 
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Figure 5. A state space portrait of plane Couette flow for Re = 400 and [L^c, Ly, Lz] = 
[27t/1.14, 2, 47t/5], projected from 61 506 dimensions to 2. The labeled points are exact equilib- 
rium (steady-state) solutions of the Navier-Stokes equation (see §[3|; the curved trajectories are 
fully-resolved time-dependent numerical integrations of Navier-Stokes projected onto the (ei, 62) 
plane defined by (|4.3|l . W^^, the Id unstable manifold of the 'lower-branch' equilibrium Ulb, 

and TzWl^u, its half-cell translation in z, are shown with thick blue lines. W^^'^\ a 2d portion 
of the unstable manifold of the 'newbie' equilibrium Unb, is shown with thin black and red spi- 
rals emanating from Unb- Similarly, the thin green lines spirally out of Uub and TzUub indicate 
and TzW^^ , the 2d unstable manifolds of Uub and its half-cell translation t^Uub within 
the S*- invariant subspace Us- Open dots along W^b show initial conditions for Newton-GMRES 
searches used to find Unb. The plane of the projection is defined in terms of the equilibrium 
solutions; it is dynamically invariant and independent of the numerical representation. See i; l4.2l 
and S I4.3l for discussions of the projection and the dynamics. 



laminar flow, one immediately, and the other after a turbulent excursion towards Uub- 
The portion of the unstable manifold of t^Ulb shown here was obtained by applying the 
Tz shift, 02 — a2, to W^g . 

We were lead to the discovery of the 'newbie' equilibrium u^b by the circular curvature 
of W^g and t^W^^ in the projection of figure El which suggested the possibihty of an 
equilibrium with a complex eigenvalue near the center of curvature. We initiated Newton- 
GMRES searches for an equilibrium at several positions between noon and three o'clock 
along Wj^B , as pictured in figure [Sj each search converged either on Ulb or on the new 
equilibrium u^b- 

The Unb unstable manifold W^^ : Within Us, the u^b equilibrium has a complex 
pair of unstable eigenvalues and one real unstable eigenvalue (figurc[4]('&j). The instability 
of the real eigenvalue is weaker than the complex pair; we omit it from consideration 
here and focus on the 2d subset W^i^'^^ corresponding to the complex pair Anb^"* with 
eigenvectors wis^'. VF^b^'^'' is shown in figure[5]as a spiral of trajectories emanating from 
Unb, calculated as discussed in 5 13.51 This simple geometric picture produces our first 
striking result: the 2d surface W^i^'^^ is apparently bounded by the Id curve W^^ . 

A heteroclinic connection from Unb to Ulb: As it approaches Ulb, VF^i^'^^ sep- 
arates along the two branches of W^^ . Since Ulb has a single unstable eigenvalue, we 
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expect that a single trajectory in W^i^''^'^ straddles the split along W^-^^ and is drawn in 
towards Ulb along its stable eigenvectors as < — *■ oo, forming a heteroclinic connection 

from Unb to Ulb. 

This is a strikingly unexpected result. In dimensions higher than two, heteroclinic 
connections are nongeneric, since it is unusual that a Id trajectory can be arranged 
to strike a particular zero-dimensional point. However, discrete symmetries and the di- 
mensi onality of the Ur,p unsta b le manifold make het e roclinic connections poss ible in this 
case (|Kevrekidis et al\1^199(j[) : iHolmes et d\ (|l996l ): ICvitanovic et all (|2008f) '). The set 



of candidate trajectories emerging from the neighborhood of u^b is increased from one 
dimension to two by the complex instability (or three if Anb is considered as well) . The 
dimensionality of state space near the target Ulb is effectively reduced to one by its 
codimension-1 set of stable eigenvalues. 

Considered in the full space U, the continuous translation symmetry increases the 
dimensionality of both the candidate trajectories and the target by two. However, the 
invariance of Us under Navier-Stokes immediately restricts possible heteroclinic connec- 
tions between the torus of u^b and Ulb translations to the four translations of Ulb within 
Ug: if a trajectory in the unstable manifold of u^b terminates at a Ulb state, it may do 
so only at Ulb, t^Ulb, t^Ulb, or Tj-^u^b- Note also that most weakly stable eigenvalues of 
Ulb, aIb* through a[b , are outside the U5 subspace, so trajectories in WCb^'^^ are forced 
to approach Ulb along the more strongly contracting eigendirections of a[b^ and a[b°^ 
(table [1]). 

The heteroclinic connection from u^b to Ulb forms a boundary between trajectories 
that decay immediately to laminar flow and those that grow towards transient turbulence. 
Those that pass near Ulb and grow to turbulence follow the unstable manifold of Ulb into 
a region near the u^b equilibrium. For Re ~ 400 and [L^, Ly, Lz] = [27t/1.14, 2, 47t/5], all 
generic initial conditions investigated so far ultimately decay to laminar. But, at higher 
Reynolds numbers and larger aspect ratios for which turbulence is sustained, wc expect 
that the u^b ^ Ulb heteroclinic connection will form a Id portion of the boundary of 
the laminar state's basin of attraction. This Id boundary should be extendable to 2c? by 
adding the third unstable eigenvalue of u^b into consideration. 

Lastly, we note that it is not possible to determine from figure [S] alone whether the 
heteroclinic connection from u^b goes to Ulb or t^^Ulb, since both of these map to the 
same point in the {61,62} plane of projection. Figure [S] (discussed below) resolves this 
question and shows that the connection is indeed from u^b to Ulb- 

Dynamics near t^Ulb and TxzUi^b- A second separation of u^b's unstable manifold 
W^NB^'^' occurs in the bottom half of figure O near TjUlb. Trajectories on the laminar 
side of t^Ulb follow its unstable manifold towards the laminar state; those on the other 
side head towards turbulence in the direction of the TjUlb unstable manifold. 

The dynamics in the region near t^Ulb can be clarified by adding a third coordinate 
03 = (u, 63) to the 2d projection of figure[5] Since 63 is antisymmetric in Xj,, the 03 coor- 
dinate distinguishes states such as TzUlb and t^^^Ulb that are related by a translation 
and so lie on top of each other in the projection of figure[51 (t^Ulb, 63) = — {txzUl,b, 63). 
Figure [H shows a 3d perspective of (ai,a2,a3) which reveals that the second separation 
of M^NB^'^', unlike the first, does not result from a heteroclinic connection between u^b 
and a translation of Ulb. As trajectories straddling the split near t^Ulb and t^j^Ulb are 
refined, they approach neither of these points. Likewise, the extensions of the refined tra- 
jectories approach neither the W"^lb = '''z^lb nor W^^lb = TxzW^^ unstable manifolds. 
The thick red trajectory in figure [S] passes closer to t^zUi^b than the other trajectories 
shown, but its recedes from t^zW^^ instead of approaching it. 
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Figure 6. The unstable manifold of Ulb and its half-cell translations W^b, t^zW^s, etc. (thick 
blue lines) and a 2d portion of the Unb unstable manifold W^g'^'^ (thin black and red lines). 
Unb is shown in only one translation (center of spiral); all four translations of Ulb are shown 
(solid dots -the unlabeled dot underneath Ulb is t^Ulb). The thick red line makes the closest 
pass to Taj^ULB of the trajectories shown. The projection is from 61,506 dimensions to 3 in the 
translation-symmetric global basis {ei, 62,63} defined by (|4.3p . 



The geometry of Wnb ' in this region is fairly complex. The separation of trajectories 
between t^Ulb and TxzUi^b suggests that another equilibrium might exist in this region; 
however, our Newton-GMRES searches initiated in this region converged on t^Ulb or 
t'xzUlb- It is clear, however, that the geometry of W^nb^'^' is shaped by the unstable 
manifolds of Ulb and two of its translations, namely, W^g, W^^^g, and W^^^^b- The 
upper-branch solution also plays a role: in figure [5] one trajectory within VFj^i^'^'' is 
drawn towards t^Uub and follows trajectories in its unstable manifold. The perspective 
of figureinialso identifies Ulb and not Ta;ULB as the endpoint of the hcteroclinic connection 
discussed above. 

Thus, with two simple state-space portraits, we have identified several regions in state 
space that trigger transitions toward qualitatively different types of flow. We expect that 
identification of such state-space regions will be extremely valuable in the development 
of nonlinear control strategies for wall-bounded turbulence. 

4.4. A local state-space portrait: the unstable manifold ofuyg in Us 

The eigenfunctions of an equilibrium provide a natural coordinate system for viewing 
its local dynamics. Within the S'-invariant subspace Ug, Uub has a single complex pair 
of unstable eigenvalues, which define a plane of local linear oscillation and two natural 
directions for a local coordinate system. The 2d portion of W^^ within Ug, which we 
denote by , departs from this plane as the distance from Uub increases and the 

magnitudes of nonlinear terms in the local Taylor expansion become nonnegligiblc. But 
since the nature of this nonlinearity was unknown, it was not immediately clear in our 
investigations how to choose a third basis function for a 3d projection of local Uub 
dynamics. We tried a variety of candidates, including principal components analysis 
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Figure 7. (a) The unstable manifolds W^^ (solid green lines) and TxW^^ (dotted black 
lines) of the upper-branch equilibrium Uub and its half-cell translation TajUuB- (b) A refined 
view of dynamics within W^^ . The coordinates {a^^a^^a^) are projections onto the basis set 
{61,62,63} that spans the plane of unstable oscillation around Uub and the direction between 
Uub and r^^UuB- See g |4.4l 



(i.e. local POD) on numerically integrated trajectories as they deviate from the plane of 
oscillation. This initial exploration suggested that the dominant nonlinear effects about 
Uub are in fact the linearized dynamics around its half-cell translation TxMvb- 

We then constructed a basis set by Gram-Schmidt orthogonalization of the plane of 
oscillation v^^ub, v^ub unstable complex eigenvalue pair Aujj^^ (see ii l3.5[) and 

(r^jUuB — Uub), that is, the direction between Uub and its half-cell translation in x. We 
indicate the Gram-Schmidt orthogonalizcd basis and coordinates with a A superscript, 
{0^,62,63} and a^{t) = (u(t),e^), to indicate its construction from the unstable Uub 
eigenfunctions and the Uub to t^Uub line. 

Figure [7] shows W^^ in the local coordinate system. In figure [7](a), W^ub^ spirals 
out from the center Uub in the plane {e^jGj} spanning v^^ubj ^')^vb.- Strong nonlincarity 
and strong trajectory separation first occur near point A: below A, trajectories continue 
the unstable linear oscillation for another cycle; above, they begin oscillation around 
TxUub, following paths similar to trajectories in TxW^^^. Figure[7|(6) shows a refinement 
of trajectories in W^ub on the upper side of the split at A. Near point B, the refined 
trajectories undergo a second split from their neighbors shown in figurc[7](a), and a third 
split among themselves at A' . This behavior and marked similarity to the trajectories 
of TxW^-i^ in figure [7](6) suggests that dynamics in this region consists of alternating 
oscillations around a symmetric pair of unstable equilibrium, in a manner reminiscent of 
the Lorenz system, until escape. 

This interpretation is reinforced by figure [H which shows the pair of W^ub'^ trajectories 
from figure [Tl'&j that split at A' together with their counterparts in t^W^^, replotted 
using the global translational-symmetric basis (|4.3p . The projection onto (ei + e2)/\/2 
and 63 was chosen because it provides a clear view of the path ABA' , and because 
these functions are symmetric and antisymmetric in Tx, respectively. Note that the two 
pairs of trajectories in W^^ and TxW^^ draw together just before A' . Each pair of nearby 
trajectories emanating from the same equilibrium splits at A' and switches allegiance with 
the pair from the opposite equilibrium, so that past A', trajectories on opposite unstable 
manifolds follow almost identical paths. The Ta;-antisymmetric long-term behavior of two 
nearby initial conditions from Uub suggests that the path from B to A' is one of weakening 
X variation, reaching small but nearly Ta;-antisymmetric x variation near A' . After A' , 
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Figure 8. A pair of nearby trajectories in the unstable manifold of Uub, from figure[7lfa^ (solid, 
green), plotted against their symmetric counterparts in the unstable manifold of r^jUuB (dotted, 
black), together with u^b and t^jUlb (blue dots, latter unlabeled), u^b and t^Uub (red dots), 
and the laminar equilibrium Ulm (black dot). Unb appears much closer to the laminar state than 
Ulb by an artifact of the projection; see figure [5] for another view. The coordinates ai, 02 are in 
the global translation-symmetric basis ei,e2 defined by (|4.3p and used in figure [5] and figure |6l 



a Tj^-antisymmetric instability comes into play, resulting in long-term ra;-antisymmetric 
dynamics. 

For the parameters of this study, the trajectories investigated so far leave the region 
of the Uub and its translations after a few oscillations, so that the Uub unstable manifold 
has the characteristics of a chaotic repeller. We expect that unstable periodic orbits can 
be found in this region, and we intend to explore this in a future publication. 

4.5. Transient turbulence 

The final stop in our stroll through plane Couette state space is an illustration of transient 
turbulence against the backdrop of the invariant structures featured in previous figures. 
For the Reynolds and cell aspect ratios studied here, all initial conditions investigated 
so far ultimately decay to laminar. Figure [9^0) shows a single trajectory, initiated as 
a perturbation of u^b, that exhibits transient turbulence and then decays to laminar 
flow. The coordinate system is ()4.3p . the same as in figure [6l The region of state space 
explored by this trajectory is typical of all observed transiently turbulent dynamics in 
Ug. The trajectory is unusually long-lived; after leaving the spiralling region around 
UfjB it wanders for some 1000 nondimensionalized L/U time units before converging on 
the laminar state, compared to more typical 200 time-unit lifetimes of other trajectories 
initialized as pertubations of u^b. 

When seen in isolation in figure [9](^aj, the turbulent trajectory shows little discern- 
able order. When plotted within the framework of invariant structures of the flow, in 
figure [9]f6j, structure is immediately evident. In this 3c? perspective, the decay to laminar 
flow is confined to a region bounded by the u^b and Ulb unstable manifolds. Transient 
turbulence occurs on the far side of laminar from these states, ai > 0.15, and in a region 
shaped roughly by the unstable manifolds of Ulb, u^b, Unb, and their half-box transla- 
tions. Close inspection shows that segments of the transient turbulent trajectory follow 
the contours of nearby unstable manifolds. 
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(b) 

Figure 9. A transiently turbulent trajectory in tiie UNsunstable manifold, (a) in isolation 
(b) in relation to Ulb, Unb, Uub, their half-cell translations, and their unstable manifolds (see 
figure [S] and figure [5] for more detailed labeling of these features) . The final decay to laminar 
of several other trajectories in the unstable manifolds of Unb and Uub are also shown. The 
projection is the same as that of figure [S] 

5. Conclusion and perspectives 

Currently a large conceptual gap separates what has been achieved for low-dimensional 
dynamical systems and the challenges we face in understanding infinite-dimensional tur- 
bulent flows. Recent computations of invariant solutions of wall-bounded flows and their 
agreement with the qualitative features of numerical simulations and experiments suggest 
that a dynamical theory of moderate- iie turbulence is within reach. We initiate a sys- 
tematic exploration of the hierarchy of exact unstable invariant solutions of fully-resolved 
Navier-Stokes equations in order to describe the spatio-temporally chaotic dynamics of 
turbulent fluid flows in terms of these states. The key advance reported here is a novel 
visualization of moderate- i?e fluid dynamics in terms of dynamically invariant, intrinsic 
and representation independent coordinate frames. The method offers an alternative visu- 
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alization of numerical and/or experimental data of any dissipative flow close to the onset 
of turbulence. In this paper, the visualizations lead to the discovery of a new equilibrium 
solution of plane Couette flow and a heteroclinic connection between two non-trivial equi- 
libria -to our knowledge the first such connection ever observed for Navier-Stokes. We 
have also computed the eigenvalues and symmetries of the three known equilibria of plane 
Couette flow in a small periodic cell with moderate Reynolds number and established 
the low-dimensionality of their unstable manifolds. 

At first glance, turbulent dynamics visualized in state space might appear hopelessly 
complex, but under a detailed examination it appears much less so than feared: it is pieced 
together from near visitations to exact coherent states interspersed by transient inter- 
ludes. Equilibria, traveling waves, and periodic solutions of plane Couette flow embody 
Hopf 's vision: a repertoire of recurrent spatio-temporal patterns explored by turbulent 
dynamics. We conceive of turbulence as a walk through a repertoire of unstable recur- 
rent patterns. As a turbulent flow evolves, every so often we catch a glimpse of a familiar 
pattern. For any finite spatial resolution, the fiow approximately follows for a finite time 
a pattern belonging to a finite alphabet of admissible fluid states, represented here by a 
set of exact coherent states. 

What new insights does the 'unstable coherent states program' offer? Normal-form 
models derived from severe truncations of spectral representations of PDEs - most fa- 
mously the Lorenz model - capture qualitatively the bifurcations and chaotic dynamics 
evocative of those observed in fluid dynamics. In contrast, exact unstable coherent states 
and periodic orbit theory should provide accurate quantitative predictions for dynamical 
observables of Navier-Stokes (such as the average turbulent drag), for a given flow, given 
flow geometry, given Re and other parameters. This description should lead to quan- 
titative predictions of transport properties of fluid flows such as bulk flow rate, mean 
wall drag, and their fluctuations. The success of computing exact eigcnfunctions and 
unstable manifolds also opens a new approach to control of turbulence in wall-bounded 
shear flows: perturbations in these directions can be used to stabilize or chaperone the 



flow t o wards a desired flui d state, and not necessarily the laminar one (jKawahara et al 
(|200,^ : IWang fit al\ (|2007[ )). 



The state-space exploration of equilibria and their global unstable manifolds presented 
here is the first step. While important in organizing the turbulent flow, equilibria, being 
static, do not actually participate in it. That role is played by the infinity of unsta- 
ble periodic orbits densely embedded in the asymptotic attractor. That it is possible to 
compute exact 3D unsta ble periodic solutions of N avier-Stokes has been de monstrated 
in the p ioneering work of iKawahara fc Kida ( 2001 ). for periodic orbits, and Viswanath 



([200761), for relative periodic orbits. However, a combination of novel and proven numeri- 
cal and analytical techniques such as variational solvers, periodic orbit theory, and group 
representation theory still needs to be developed in order to systematically explore the 
hierarc hy of such solutions and to derive the statistics of the flow through periodic orbit 
theory (jCvitanovic et all (|2007t )'). 
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for support. J.H. thanks R. Mainieri and T. Brown, Institute for Physical Sciences, for 
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Table 1. (left) Least stable eigenvalues of the laminar equilibrium Ulm for 
[Z/x, I/y, -Lz] = [27t/1.14, 2, 47t/5] and Re = 400, computed by Arnoldi iteration, compared to 
Stokes (S) and heat-equation (H) eigenvalues from analytic formulas. This serves as a test 
of accuracy for our channelflow.org codes. The heat-equation eigenfunctions have the form 



u(x,t) 



' sm{nkyy/2) cos{2Ttkzz/ L~)) a for ky even and e cos(7tA:yy/2) cos(27tfcz2;/Z/z)) x 



for ky odd, and eigenvalues A — —{n ky/4:-\-ATrk'j,/L^)/Re. The Stokes eigenvalue listed is 
the low e st-ord er mode with v component of the form v(y)cos{2nz/ Lz), v{y) even in y (see 
IWaleffd (|1997l ')'). The eigenvalues are ordered in the table by decreasing real part. All laminar 
eigenvalues are real, (right) Ulb equilibrium stability eigenvalues A = fi ± iuj and symmetries 
of corresponding eigenvectors at same parameter values. The zero eigenvalues result from the 
continuous translation symmetry of the flow. 
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Table 2. Stability eigenvalues X — /iiiz ico and symmetries of corresponding eigenvectors: 
(left) Unb, (right) UuB equilibrium for [Lx,Ly,Lz] — [27t/1.14, 2, 47t/5] and 7?e — 400. 



Appendix A. Tabulation of numerical results 

Tables [Hand [2] list the least stable linear stability eigenvalues of the Ulm, Ulb, u^b, and 
UuB equilibria, together with symmetries of corresponding eigenfunctions. The unstable 
eigenvalues together with a set of the least contracting stable eigenvalues are also shown 
in Figure [4l All numerical results tabulated in this appendix are computed for plane 
Couette flow with Re = 400 and [La;,Ly,L^] = [27t/1.14, 2, 47t/5]. FuU sets of exact 
invariant solutions, their linear st ability eigenvalues and eigenfunctions arc available on 

or can be obtained by a request to authors. 



channelf low . org, 
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